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r^ . Abstract 

We present results of a Monte Carlo simulation study of the hadron energy 
response (for 1 GeV < E < 15 GeV) of the magnetized iron calorimeter detector, 
^ . ICAL, proposed to be located at the India-based Neutrino Observatory. Using a 

^'* \ GEANT4 modeling of the ICAL, interactions of atmospheric neutrinos with tar- 

get nuclei are recorded. We discuss a method of calibration of the hadron energy 
10 ■ using the hadron hit multiplicity in the active detector element. We first study 

'!::j- \ the detector response with single pions propagating through the detector. Then 

C3 ' the average response of hadrons produced in atmospheric neutrino interactions is 

analyzed using NUANCE-generated neutrino events. The shape of the distribution 
and an appropriate fitting function are examined in detail. Finally, the hadron en- 
ergy resolution is determined as a function of energy. We find an energy resolution 
r> ■ of around (60-40)% for hadron energies in the range 2-15 GeV. 
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1 Introduction 

The India-based Neutrino Observatory (INO) [1] is a planned underground laboratory in 
South India with the long term goal of conducting experiments on cosmic rays, neutrinos 
and astrophysics. The primary focus of INO in the first phase is the study of atmospheric 
neutrinos with a magnetized iron calorimeter (ICAL) detector similar in concept to the 
design of the MONOLITH detector [2]. The detector is designed to observe neutrino (and 
anti-neutrino) interactions in the GeV range. The main aim of the ICAL detector is a 
precise measurement of the oscillation parameters both for neutrinos and anti-neutrinos; 
recent measurements of the mixing angle 613 in reactor experiments [3] |1] [S] will enable 
the ICAL to pin down the mass ordering of neutrinos by separate measurements of i>n 
and z/^ interactions by exploiting matter effects in the earth. 

The specifications of the magnetized ICAL detector, consisting of three identical mod- 
ules with layers of iron plates interspersed with resistive plate chamber (RPC) detectors 
and approximately 50 kton in mass, are given in Table [H RPC detectors have been 
chosen as the active detector elements by virtue of their good position (cm) and timing 
(ns) resolution, high detection efficiency for charged particles, long term stability, and 
coverage of substantial area. The RPC readout strip width for the present simulation 
study has been assumed as 1.96 cm. ICAL is expected to have good charge identifica- 
tion efficiency and good tracking and energy resolution, especially for muons produced 
in charged current (CC) interactions of the detector material with atmospheric muon 
neutrinos. A schematic illustration of the detector is shown in Fig. [1] 



No. of modules 




3 


Module dimension 


16 m 


X 16 m X 14.5 m 


No. of RPC layers 




150 


Iron plate thickness 




5.6 cm 


RPC unit dimension 


1.84 m 


X 1.84 mx 24 mm 


Total no. of RPC units 




28800 



Table 1: ICAL detector specifications. 



Atmospheric neutrinos (antineutrinos) interact with the iron target in the detector 
through quasi-elastic (QE), resonance (RS) and deep inelastic scattering (DIS) processes 
as shown in Table [2], where / is the charged lepton, A^, A^' are nucleons, while X denotes 
any set of possible final hadrons, with typically one pion and a nucleon in RS, and 
multiple hadrons in DIS interactions. At these energies, coherent interactions and the 
interactions with the electrons in the detector are rare and are ignored. 

In this paper, we are concerned mainly with the charged current interactions in 
which, apart from a muon, there will be hadrons in the final state. While the final 
state muon in the charged-current z/^ interaction is reconstructed well for its energy and 
direction, the hadrons in the final state are observed mainly as a bunch of hits close 
to the interaction vertex. The physics analysis crucially depends on the precision with 
which the neutrino energy is reconstructed, and hence on the precision of reconstruction 
of the energies of the hadrons and muon. The hadrons consist mainly of pions (about 




Figure 1: Left: Schematic of ICAL detector. Right: a close up view of the stack and 
RFC insertion method is shown. 



Frocesses 


cc 


NC 


QE 
RS & DIS 


ui N ^l- N' 
{-) <+^ 


^u} N -^ ViN' 
ui N -)■ viX 



Table 2: The neutrino interaction processes. 



85% of events on the average), kaons, and also nucleons, including the recoil nucleon 
which cannot be distinguished from the remaining hadronic final state. Both neutral 
and charged pions may be produced. The neutral pion decays immediately, giving rise 
to two photons, while the charged pions propagate and develop into a cascade due to 
strong interactions. The main uncertainty in determining the incident neutrino energy 
comes from the uncertainty in estimating the energy of these hadrons. Typically a muon 
leaves a distinct track whereas the hadrons leave a shower of hits close to the vertex. 
A visualization of a neutrino DIS event with large hadron energy component in the 
simulated ICAL detector (using the VICE event display package [B]) is shown in Fig. [21 
For the interaction v^N — )■ /iX, the incident neutrino energy is given by 



E„ 



E„ + El 



had' 



(1) 



where 'El^^^ is the difference in the total energy of the final state hadrons (including the 
recoil nucleon) and the energy of the initial state nucleon — including the rest mass, 
Fermi momentum and binding energy. We use Ej^^^ as an observable for hadron energy 
since it would enable us to reconstruct the neutrino energy in a straightforward manner. 
The analysis is done in two parts. Mono-energetic charged pions are propagated 
through the detector and their characteristics studied. These are used to determine 
a suitable probability distribution function (FDF) at each energy. Next, we consider 
hadrons produced through neutrino interactions. Here in principle all types of hadrons 
contribute to the energy, though it is dominated by pions at the energies of a few GeV. 
E^^j can be reconstructed in ICAL by devising a calibration with respect to the number 



of hits in a hadron shower. In general the energy resolution of hadrons is affected by 
various factors like shower energy fluctuation, leakage of energy and invisible energy loss 
mechanisms. 
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Figure 2: A deep inelastic neutrino interaction event in the simulated ICAL detector 
using VICE event display package. 

In this paper we use a full simulation of the event propagation in ICAL to assess 
the detector's capability to estimate hadron energy. In section 2, we first discuss the hit 
distributions of hadrons in the detector. This is done by a Monte-Carlo (MC) generator 
that generates hadrons at fixed energies up to 15 GeV. The hits are are estimated using 
the hit multiplicity information of events generated with appropriate position- and angle- 
smearing. We then specifically discuss the pion energy (E^) resolution by propagating 
pions with fixed energies through the simulated detector, focusing on the appropriate pdf 
to be used for analyzing the pion hits. In section 3, we use the hadrons (not just pions) 
generated in atmospheric neutrino-nucleus interactions using a neutrino event generator 
NUANCE pj, and study the net hadron energy resolution as a function of Ej^^^^. In 
section 4, we discuss the calibration of Ej^^^^ using the hit information. We conclude in 
section 5. 



2 Energy resolution with fixed energy pions 

The ICAL detector, with the parameters specified in Table [1], is simulated using the 
GEANT4 package [TJ [8]. For more details about the coding of the detector geometry 
and the nature of hits generated, see Ref. [10]. Our focus in this section will be on the 
analysis of fixed-energy single-pion hit distributions. 

When a charged particle propagates through the ICAL detector, hits in the X and Y 
strips of the RFC layers are recorded [12]. The layer number provides the z-coordinate. 
Thus the full position information is available. A muon usually leaves one or two hits 
per layer and so the hits from both X and Y strips can easily be combined to obtain the 



number of hits and their position coordinates (x, y) in a given layer. However, in the 
case of hadron shower there are multiple hits per layer, and combining X and Y strip hits 
leads to some false count of hits (ghost hits). To avoid the ghost hit counts, the energy 
calibration may be done with counts from either X or Y strips. The variables x-hits 
and y-hits store the number of hits in the X and Y strips of the RPC, respectively. The 
maximum of x-hits or y-hits is stored as the variable orig-hits. In Fig. |31 the comparison 
of these three types of hit variables are shown for tt^ of energy 3 GeV. As is clear from 
Fig. [3l we could have used any of the x-, y-or orig-hits in the analysis; we choose to use 
orig-hits as the unbiased parameter. 
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Figure 3: The comparison of the distributions of x-hits, y-hits and orig-hits for (left) it 
and (right) 7r+ of energy 3 GeV. 

Unless otherwise specified in this section, the total number of events for each input 
energy value is 10000. Each event is randomly generated over a volume 200 cm x 200 cm 
X 200 cm in the central region of the ICAL detector. In addition, since there is very little 
impact of the magnetic field on the showers produced by hadrons, the hadron direction is 
uniformly smeared over zenith angle < ^ < tt and azimuth of < < 2it. The angles 
are denoted with respect to a reference frame, where the origin is taken to be the center 
of the detector, the 2;-axis points vertically up, while the plates are horizontal in the x-y 
plane. This serves to smear out any angle-dependent bias in the energy resolution of the 
detector by virtue of its geometry which makes it the least (most) sensitive to particles 
propagating in the horizontal (vertical) direction. 

In Fig. mthe hit distributions in the detector for pions, kaons, and protons at various 
energies in the range of 1 to 15 GeV are shown. We observe that for all these hadrons the 
hit patterns are similar, though the peak position and spread are somewhat dependent 
on the particle. The detector cannot distinguish the specific hadron that has generated 
the shower. Since hadrons produced in neutrino interactions with ICAL are primarily 
charged pions, we focus our attention on the detector response to pions in this section 
and consider a more general admixture of different hadrons in Section 3. However, since 
the energy response of all hadrons appears to be rather similar, there should be no 
significant information loss even when we cannot identify the final state hadrons. 

We will therefore use the hit distribution of charged pions at various energies to 



determine the detector response. In particular we determine the mean and variance 
of the distribution and their energy dependence. Later on we will discuss a suitable 
probability distribution function for the hit patterns generated. This information is 
required in order to determine the energy calibration for use in physics studies with 
ICAL. 
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Figure 4: The hit distributions at various energies (angle averaged) for tt^, 7r°, K"^ and 
protons propagated from vertices smeared over the chosen volume. 



2.1 Analysis of the pion (tt ) hit pattern 

The pion hit distributions at sample values of E = 3, 8 GeV are shown in Fig. 0, where 
the arithmetic mean and rms value obtained from the histogram are also listed. Typical 
patterns show a mean of roughly 2 hits per GeV as seen from the figure, but with long 
tails, so the distribution is not symmetric. In addition, several events yield zero hits in 
the detector at lower energies, such events are virtually absent at higher energies. 

Fig. [5] also shows the fits to the Gaussian distribution. As can be seen from the figure, 
the quality of the fit at the lower energy (E = 3 GeV) is poor, with a y^ per degree 
of freedom, x^/ndf ~ 4.6. This improves to y^ /xidi ~ 2.2 for E = 8 GeV. This clearly 
indicates that the hit pattern is non-Gaussian for lower energy pions. This is also easily 
understood because of the long non-Gaussian tails. However, at high energies the hit 
distribution is well approximated by a Gaussian. 



6 



Mean. 




3GeV 



6.019 ± 0.035 
3.484 ± 0.025 



Mean, 



X^ndf =118.6/26 



Meanr- 



= 6.022 ± 0.043 
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Figure 5: The hit distributions at 3 GeV (left) and 8 GeV (right), for pions propagated 
from randomized vertices over a volume of 200 cm x 200 cm x 200 cm in the detector. 
A Gaussian fit to the histogram is included for later convenience. 
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Figure 6: Mean hit distribution as a function of pion energy with 5.6 cm iron plates. Left 
(right) panel shows the mean hit for t\~ (tt"^). Black dots are the results of simulation 
while the curves are fits to this simulated data. 



The arithmetic mean n(E) of the number of hits increases with increasing pion energy, 
and saturates at higher energies. It may be approximated by 



n(E) = no[l-exp(-E/Eo)], 



(2) 



where uq and Eq are constants. This fit has to be used with some care, since uq and Eq 
are sensitive to the energy ranges of the fit. In Fig. [6], we compare the means of the hit 
distributions for tt" and vr"*" events in the detector, in the energy range 1-15 GeV. We 



also show the fits to Eq. 



As expected, the mean number of hits for tt and vr"*" are 



identical within expected statistical fluctuations. The value of Eq is found to be ~ 30 
GeV. Since the energies of interest for atmospheric neutrinos are much less than Eq^ we 
may use Eq. ([2]) in its approximate linear form ii(E) = UqE/Eo. Fig. Elalso shows a fit 
to this linear form. 



In the linear regime (E <^ Eq), we have 

n(E) 



no 



E 
Eo' 



(3) 



and therefore the energy resolution may be written as 

An{E) 



a 
E 



n(E) ' 



(4) 



where An is the RMS width of the distribution. In the rest of the paper we shall use the 
notation a/E for energy resolution, and use Eq. (jlj) for the rest of the analysis. 
The energy resolution of pions may be parametrized by 



a 

E 




(5) 



where a and b are constants. In the literature a is often referred to as the stochastic 
coefficient which incorporates both statistical and systematic uncertainties. In Fig. [TJ 
we show the vr" and vr^ energy resolutions as functions of pion energy, and extract the 
parameters a and 6 by a fit to Eq. (|5]) over the energy range 1-15 GeV. It is seen that 
while the fits to both vr" and 7r+ yield virtually identical values for a and b, the quality 
of fits is rather poor at low energies. 
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Figure 7: Energy resolutions for ir (left) and 7r+ (right) for different energies shown as 
a function of pion energy. Also shown is the fit to this resolution according to Eq. ([5]). 



As observed in this section, the Gaussian distribution is not a good approximation 
at low energies, while, it approximates the hit distribution well at higher energies. It is 
desirable to find a suitable fit function that correctly reproduces the mean and variance 
at all energies, and is well-approximated by a Gaussian distribution at higher energies. 
We shall do so in the next section. 



2.2 A suitable probability distribution function 

From Fig. HJ it is clear that the hit distribution approaches the Gaussian pdf only at 
energies greater than around 5 GeV. Since the bulk of atmospheric neutrino events, the 
ones that ICAL detector will record, are concentrated in the energy range less than about 
10 GeV, it becomes important to find a suitable pdf for the hit distribution. This should 
necessarily agree with the Gaussian pdf at high energies. We consider this question in 
detail in this section. 

Here we analyze the hit distributions using three fitting functions. As we have seen, 
the hit distribution is asymmetric, particularly at lower energies. We therefore compare 
the following functional forms: (i) Gaussian, (ii) Landau convoluted with a Gaussian, 
and (iii) Vavilov [TT]. The energy loss of a particle passing through moderately thick 
absorbers is commonly modelled using the above three distributions depending on the 
energy of interest. 

In Fig. [HI we show the sample hit distributions for 10000 pions propagated through 
the detector from randomised vertex points, at sample energies of 3 GeV and 8 GeV, 
fitted to the three functional forms. A comparison of the fit parameters is shown in 
Table [3], where we also indicate the P-values obtained from the Kolmogorov-Smirnov 
(KS) test. It is clear that the best fit is obtained with the Vavilov pdf. A sweep of all 
other energies shows that the Vavilov pdf consistently gives a better fit than the other 
distributions. This is what we will use in the analysis that follows from now on. We 
however note that even though the x^ differs substantially, the mean and a obtained 
from different fits are quite similar. We shall quantify this in section 12.41 



E. 


Parameters 


Histogram 


Gaussian 


Landau 


Vavilov 


3 GeV 


xVndf 

P (KS test) 

n 


6.019 


4.6 

8.1 X 10-6 

6.022±0.043 


4.2 

7.6 X 10-6 

6.021±0.043 


1.7 

0.95 

6.018±0.044 




a 


3.484 


3.671±0.042 


3.601±0.048 


3.431±0.052 


8 GeV 


xVndf 

P (KS test) 

n 


14.813 


2.204 

2.6 X 10"^ 

15.151±0.061 


2.01 

1.2 X 10"^ 

15.062±0.063 


1.607 

0.92 

14.781±0.063 




a 


6.114 


6.138±0.051 


6.129±0.053 


6.124±0.068 



Table 3: Comparison of the three different pdf fits for 3 GeV and 8 GeV pions. The 
X^ values as well as the Kolmogorov-Smirnov (KS) statistics quantify the substantial 
differences in the quality of fits with different pdfs. Note that the mean and sigma are 
less sensitive to the choice of pdf. For comparison the mean and rms obtained directly 
from the hit distribution are also shown in the third column. 
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Figure 8: The hit distributions at sample energies 3 GeV and 8 GeV, shown with statis- 
tical errors and fitted with (top) Gaussian, (middle) Landau convoluted with Gaussian, 
and (bottom) Vavilov distributions. The resulting fit parameters and quality of fits are 
shown in Table [31 



10 



2.3 Energy resolution using the Vavilov fit 

We now parameterize the energy resolution of pions using the Vavilov fit, which has four 
parameters, Pq, Pi, P2 and P3 (see Appendix [K\ for details). In Fig. |9l we show the 
variation of the Vavilov parameters Pq, Pi, P2, P3 with pion energy. The mean and a are 
obtained from these parameters via 



Mean = (7 - 1 - InPo - Pi) P3 + P2; a^ = (2 - Pi) 



2Pr 



(6) 



where 7 is the Euler's constant. The main features of the Vavilov distribution may be 
captured by the two parameters mean and sigma. However the remaining information 
in the Vavilov parameters enables a good fit even in the regions of low hit multiplicity 
and the tail of high hit multiplicity. 
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Figure 9: Variation of the Vavilov fit parameters Pq, Pi, P2, and P3 with pion energy. 



Notice that beyond E ~ 6 GeV, we have Pq 
to a Gaussian. 



10, so that the Vavilov distribution reduces 



Note that for Pq > 10, Vavilov distribution reduces to a Gaussian pdf with the 
corresponding mean and a. In this analysis this happens for energies greater than ~ 6 
GeV. The actual power of Vavilov distribution is therefore apparent in the energy range 
1-6 GeV, which is the main region of interest for hadrons in the ICAL detector. 

Thus having extracted the mean and a as discussed above, in Fig. [10] we plot the 
resolution a/E as a function of energy. This is then fitted to Eq. ([5]) to obtain the 
coefficients a and b. The fit has been done separately for the energy ranges 1-5 GeV, 
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Figure 10: 
resolution 



Resolution cr/E plotted using the Vavilov fit for as a function of energy. The 
is fitted to Eq. ([5]) in three different energy ranges as shown. 



5-15 GeV, and the combined energy range 1-15 GeV. The low energy region (1-5 GeV) 
is the region where Vavilov pdf gives the best fit, while the Gaussian distribution fails. 

2.4 Comparison of the means and cr's for different pdfs 

Before we conclude the analysis of fixed energy single pion events, a few remarks are in 
order about the various fits used. Though the Vavilov distribution gives a good fit to 
the data, it has four correlated parameters, as opposed to the two parameters meancauss 
and (JGauss in the Gaussian pdf. To a first approximation the meanvaviiov and avaviiov are 
sufficient to determine the energy resolution. Often, for some applications, for example 
where detector parameters are checked for performance, a simpler alternative may be to 
use the fit-independent parameters meauArith and aArith, the arithmetic mean and the rms 
width, respectively. In Fig. [HI the mean and a obtained directly from the hit distribution 
is compared with those obtained from the Gaussian and Vavilov distributions. 

The plots in Fig. [11] show that meauArith and the means obtained from the dis- 
tributions are comparable with marginal differences. However, the comparison of a's 
shows that at low energies acauss deviates from (and is higher than) (XArith- This can be 
clearly seen in Fig. [121 where the ratios of the fitted mean and a with respect to the 
corresponding fit-independent values are plotted as functions of energy. The marginally 
higher value of meangauss as compared to meauArith at higher energies is the result of 
the long tail in the hit distribution histogram. Also acauss exceeds the others at lower 
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Figure 11: The comparison of mean (left) and a (right) obtained from the Vavilov and 
Gaussian pdfs, as well as the fit-independent quantities. 
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Figure 12: The ratio of mean and a from Vavilov and Gaussian distributions to meauArith 
and ciArith- 

energies because of the many events with zero hits at lower energies. Interestingly, at 
higher energies, while acauss and avaviiov agree, they marginally underestimate the UArith- 
This is again due to the presence of the tail. While care is required at low energies, the 
resolution is fairly distribution independent for higher energies, beyond 6 GeV. 

3 Energy resolution with NUANCE-generated neu- 
trino events 

The previous section contained an analysis of the energy resolution with single pion 
events. In actual practice we are interested in the energy resolution of hadrons produced 
through atmospheric neutrino interactions. In this section, we analyze the energy reso- 
lution of hadrons produced in charged-current (CC) atmospheric z/^ interactions in the 
detector. Both atmospheric neutrino {u^) and anti-neutrino (l?^) events in ICAL are 
generated using the neutrino event generator NUANCE (v3.5) [9]. Though these are 
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mainly pions, there are some events with kaons (about 10%) and other hadrons. There 
are in general multiple hadrons, especially in the DIS events. As discussed earlier, it is 
not possible to distinguish among different hadrons from the "hadron shower" seen in 
the hit pattern. However, since the hit distribution of various hadrons are similar to each 
other (see Fig. H]), and the NUANCE generator is expected to produce an approximately 
correct mixture of different hadrons at all energies, it is enough to determine the hadron 
energy resolution at ICAL through an average of NUANCE events, without having to 
identify the hadrons individually. 

The analysis is done with 1000 kt-years of data (equivalent to 20 years of exposure 
with the 50 kton ICAL module). The events are angle-averaged and binned in the hadron 
energy E'-^^^^ = E,^ — E^. The hit distributions (of the variable orig-hits) are analyzed as 
before, in each energy bin, averaged over all angles. We show the meauArith and o"Arith of 
these distributions in Fig. [131 As expected, these are similar to the mean and a obtained 
earlier with fixed energy pions. Note that while the energy is kept fixed to a single 
value in the latter case, the events are binned in energy in the NUANCE-based analysis. 
Since the mean hits grow approximately linearly with energy, we use the same linearized 
approximation as before to obtain the energy resolution cr/E = An/n. 
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Figure 13: MeauArith and (XArith of the hit distribution as a function of hadron energy 
Ejj^jj. The bin widths are indicated by the horizontal error bars. 



3.1 Analysis of the hit distribution 

Now we use the Vavilov distribution function as a fit to the hit distribution of hadrons. 
The hit distributions and their Vavilov pdf fits for selected energy ranges are shown in 
Fig. [m While performing the fit, we use the central value of Ej^g^^. 

In Fig. [15] we show the energy dependence of the parameters Pq, Pi, P2, and P3 of 
the Vavilov distribution (see Appendix [X] for details on these parameters) obtained from 
the fits to the NUANCE data. It can be seen that these results are similar to the fits 
obtained from the single pion analysis. The changes, mainly a small increase in Pq, P2 
and P3 at energies < 6 GeV, are primarily because the hadrons in a bin do not have 
identical energies like the fixed-energy pions in the earlier analysis did. 
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Figure 14: Hadron hit distributions in different E^,^^ energy bins fitted to the Vavilov 
distribution function in each energy bin. 

The four parameters P, are sufficient to specify the distribution completely and may 
be directly used in the physics analysis. However, the quantities meanvaviiov and uvaviiov 
are more useful in intuitively assessing the performance of the detector. The dependence 
of the meanvaviiov and avaviiov on energy is shown in Fig. [161 
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Figure 15: Variation of the Vavilov fit parameters Pq, Pi, P2, and P3 with E^^^^ from 
NUANCE data. The horizontal errors indicate the bin widths. 
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Figure 16: The variations of meanvaviiov and avaviiov with hadron energy for NUANCE- 
generated data. 
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3.2 Dependence of the energy resolution on neutrino oscilla- 
tions 

Note that the effective energy response we obtained from the NUANCE-generated data 
is an average over the mixture of many hadrons that contribute to hadron shower at 
all energies. The fractional weights of different kinds of hadrons produced in neutrino 
interactions may, in principle, depend upon neutrino oscillations. In this section we 
check if there is any bias in the hadron energy resolution due to oscillations. We use the 
standard oscillation parameter best-fit values [13] for the mixing angles and mass-squared 
differences given in Table H] to generate the events with oscillations. The resolutions 
obtained without and with oscillations are compared in Fig. [T71 As expected, the 
resolutions are virtually insensitive to oscillations. The marginally higher values of a/E 
with oscillations that is seen in the figure is an artifact of the smaller number of events 
with oscillations. 
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Table 4: The values of oscillation parameters used in the analysis. 
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Figure 17: Comparison of energy resolution a/E of hadrons from charged-current muon 
neutrino events from oscillated and unoscillated events. 



4 Hadron energy calibration from the number of hits 

In the previous sections, we have shown the dependence of hadron energy resolution in 
ICAL on hadron energy. When actual data become available, the hit multiplicity will 
be known, and the hadron energy E'^^^ will have to be calibrated. We therefore calibrate 
^had with respect to the hadron hit multiplicity in a simulated event. 
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The calibration of the hit muhiphcity in terms of the hadron energy not only depends 
on the detector response discussed above, but also on the energy distribution of hadrons 
in the detector (that is, the sample of events from oscillated neutrino interactions). To 
take care of this we have used atmospheric neutrino events generated with the best fit 
values of the oscillation parameters, as in Table IH as input to the simulation code. We 
have checked that changing the oscillation parameters within their current la allowed 
ranges does not change the results appreciably. 

The simulated data were divided into hit_n bins, where hit_n corresponds to n hadron 
hits. Note that the same number of hits can arise from events with different hadron 
energies. The distribution of these energies was obtained for each bin, and fitted with 
the Vavilov distribution. In Fig. [TSl the distribution of energy for the bins with 10 and 
15 hits are shown. 
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Figure 18: The distribution of hadron energy for (left) 10 hits, and (right) 15 hits. 

Fig. [19] shows the calibration plot, with the mean E^^^^ and the error o obtained from 
the Vavilov fit. Now from the hit information of the hadron shower in any event in ICAL, 
we can estimate the hadron energy using this calibration plot. This can then be used to 
reconstruct the energy and direction of the incident neutrinos. 
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Figure 19: Calibration plot of Ej^^^^, where mean and o from the Vavilov fits are repre- 
sented by the dots and error bars, respectively. 



5 Summary and Conclusions 

The ICAL detector is expected to be a good tracking detector for the muons produced 
in charged current atmospheric u^ events. However, in order to reconstruct the neutrino 
energy, the estimation of the energy of the hadrons produced in such interactions is 
crucial. Using a GEANT4-based simulations framework of the ICAL detector, we have 
explored the energy calibration of hadrons from the hadron hit pattern and obtained the 
hadron energy resolutions. 

We now summarize our main results. We have found that the Vavilov distribution 
function is a good functional form to fit the hadron hit distribution. Analyses, first with 
fixed-energy pions, and later with a mixture of hadrons from atmospheric Ufj_ interaction 
events, show that a hadron energy resolution in the range (60-40)% at 2-15 GeV is 
obtainable for hadrons produced in charged- current interactions. We have presented the 
hadron energy resolution as a function of the incident hadron energy, as well as the 
calibration plot for the hadron energy, given the hit multiplicity in the detector. 

In principle, it is also interesting to study the angular dependence of these resolutions. 
This will enable the direction of the shower to be determined and will add more infor- 
mation to the reconstruction of the neutrino kinematics. The methodology followed and 
the issue of the reconstruction of the hadron shower direction will be discussed elsewhere 
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A The Vavilov distribution function 

The Vavilov probability density function in the standard form is defined by 

1 /"C+ioo 

P(x;n,P') = — (^{s)e^^ds, (A.l) 

^VrZ J c—ioo 

where 

(j){s) = e^e^^'\ C = k{1 + (3^-i) (A.2) 

and 



iIj{s) = slnfi; + (s + /3^k) ■ 



1 

dt - 7 

t ' 

,0 
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-Ke"^/" (A.3) 



where 7 = 0.5772156649 ... is the Euler's constant. The parameters mean and variance 
{a'^) of the distribution in Eq. (]A.1|) are given by 



2- (3^ 
mean = 7 — 1 — In k — /3^; cr^ = . (A. 4) 

For K < 0.05, the Landau distribution approximates Vavilov very well while for n > 10, 
the Gaussian approximation works well. 

We have used the Vavilov pdf which is built into ROOT with the mean shifted from 
zero. In order to fit the desired distribution the following redefinition of parameters 
has been adopted: The original Vavilov pdf, that is, P(x; n, /3^) has been modified to 
(P4/P3) P(x — P2/P3; Pq, Pi)- These changes account for normalization (P4) and the 
shift of the mean to a non-zero value. This gives the modified mean and variance as 

p2 

Mean = (7 - 1 - InPo - Pi) P3 + P2; a'' = {2 - F^) -f . (A.5) 

The errors on mean and a therefore arise from highly correlated parameters Pj. These 
have been determined from standard analytical formulae |15] . 
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